Time-dependent Ginzburg–Landau equations for multi-gap superconductors
Li Minsi, Gu Jiahong, Du Long, Zhong Hongwei, Zhou Lijuan, Chen Qinghua
School of Science, Guangxi University of Science and Technology, Liuzhou 545006, China

 

† Corresponding author. E-mail: chenqinghua2002@163.com

Project supported by the National Natural Science Foundation of China (Grant Nos. 11564003 and 11865005) and the Natural Science Foundation of Guangxi Province of China (Grant No. 2018GXNSFAA281024).

Abstract

We numerically solve the time-dependent Ginzburg–Landau equations for two-gap superconductors using the finiteelement technique. The real-time simulation shows that at low magnetic field, the vortices in small-size samples tend to form clusters or other disorder structures. When the sample size is large, stripes appear in the pattern. These results are in good agreement with the previous experimental observations of the intriguing anomalous vortex pattern, providing a reliable theoretical basis for the future applications of multi-gap superconductors.

1. Introduction

Conventional superconductors can be divided into two types (type-I and type-II) according to their different responses to an external magnetic field. Both the types of superconductors expel weak magnetic fields. In bulk type-I superconductors, the superconductivity abruptly disappears when the magnetic field exceeds a certain critical value. However, in type-II superconductors, the superconductivity and magnetism can coexist when the magnetic field is between two critical values. In this case, the magnetic field penetrates into the material, resulting in the formation of microscopic magnetic vortices, and the regions outside these microscopic vortices will remain in the superconducting state. This kind of coexistence also occurs in type-I superconducting thin films when the magnetic field is strong enough but smaller than the critical field.[1,2] In the Ginzburg–Landau (G–L) theory, these two types of superconductors are differentiated by the G–L parameter, κ = λ/ξ, where λ is the penetration depth of the magnetic field and ξ is the coherence length.[1,3,4] The G–L parameters also determine the interaction between the vortices.[5] In type-I superconductors, the G–L parameters are smaller than . When the type-I superconducting sample is a thin film, the vortices generated by the magnetic field will attract each other and form the macroscopic normal magnetic domains.[2] On the contrary, for type-II superconductors, the G–L parameters are larger than , and the vortices repel each other, leading to the presence of vortex lattice. If the G–L parameters are exactly equal to , there is no interaction between the vortices.

Recently, a novel type of superconductors called multi-gap superconductors have attracted a great deal of attention from researchers.[69] In the multi-gap superconductors, the interaction between vortices becomes more complicated. The multi-gap superconducting states are characterized by multiple coherence length, ξ1 , ξ2, …, ξN, and only one penetration depth λ. It is possible that . In other words, type-I and type-II superconducting states coexist in the same sample. It is proposed theoretically that the vortices may repel each other at short distance and attract in long range. This nonmonotonic inter-vortex interaction could lead to novel magnetic responses in the superconductors.[1013] The novel magnetic responses, including unconventional disordered magnetic vortex patterns, cluster-like structures in small-size samples and stripes/gossamers in larger ones, were firstly experimentally observed in a two-gap superconductor, MgB2,[1416] and then in Sr2RuO4, LaPt3Si and many other materials.[1726] Moshchalkov et al termed this phenomenon as type-1.5 superconductivity.[14] Meanwhile, extensive theoretical investigations have been carried out on this novel phenomenon.[6,7,1013,27,28] However, to the best of our knowledge, no theoretical results are in good agreement with the experimental data, especially in vortex scale. This means that it is necessary to find a new method to calculate the novel vortex pattern in such superconducting systems. It is lucky that another important method, time-dependent Ginzburg–Landau theory (TDGL), which is normally used to investigate the dynamics of the vortices in superconductors, has been wildly accepted to understand the formation mechanism and time evolution of vortex states in single-gap samples.[2934] This is because the TDGL is consistent with an equation of the order parameter similar to the Schrödinger equation in quantum mechanics, and can be derived from the microscopic BCS theory, and thus it is more reliable than other simulation methods.

In this paper, we aim at providing a theoretical understanding of the formation of above-mentioned unconventional vortex patterns by solving TDGL equations using finite-element techniques for multi-gap superconductors. Our real-time simulation of the formation and evolution of magnetic vortex patterns is in very good agreement with the previous experimental observations,[1416] which has not been reported before. The theoretical method as well as simulation results presented here may have important implications for the future design of superconductor vortex-based devices. The paper is organized as follows. In Section 2, we present the details of our numerical formalism. In Section 3, the simulation results are discussed and compared to the experimental ones. All findings are finally concluded in Section 4.

2. Theoretical model

Unlike a single-gap superconducting system, in two-gap situation, one must consider two order parameters to describe the superconducting properties. One of the first Ginzburg–Landau descriptions of multi-gap superconductors was developed by Zhitomirsky and Dao, starting from microscopic theory.[35] Combining with the consideration of the Josephson coupling between the two gaps, resulting from the tunneling of the Cooper pairs from one gap to the other, the free Ginzburg–Landau energy functional is expressed as[35,36]

where αi and βi are the phenomenological Ginzburg–Landau coefficients, and i = 1,2 is the gap index, , and e*, respectively, are the effective mass and charge of the Cooper pairs, c is velocity of light in vacuum, A is the vector potential associated to the external applied field, ψi are the order parameters of the two gaps, and γ is the Josephson coupling strength. For MgB2, γ > 0, while for iron-based superconductors, γ is presumably negative.[37] To make the equations dimensionless and then to minimize the F with respect to the ψ1, ψ2, and A, we obtain the GL equations as follows:[28]

The single-gap TDGL equations and their solutions have been well discussed in Refs. [29,30]. The equations are

where D is the diffusion constant, σ is the normal conductivity, ϕ is the applied electric potential, and Js is the supercurrent density with . The dimensionless forms are as follows:

Similarly, for the two-gap TDGL, we can also have the dimensionless forms as

It is obvious that equations (9)–(11) are difficult to solve directly. In order to integrate them, we generalize the numerical method developed for single-component superconductors[29,30] to two-component superconductors and solve the TDGL equations numerically. We use the link variable defined as[2830]

with μ = x, y. Here the magnetic field H is in the z direction. In order to use the finite-element technique we must first discrete the system. Figure 1 is the configuration of the lattice used in simulations. The lattice points are denoted by i, j, k, … and the lattice constants in x and y directions are the ax and ay, respectively. Suppose ax = ay, then the TDGL equations can be rewritten as

where and ;

where and .

Fig. 1. Configuration of the lattice used in simulations. The lattice points are denoted by i, j, k,… and the lattice constants in x and y directions are the ax and ay, respectively.
3. Simulations and result discussion

In this work we configurate the sample grid as 64 × 64, 128 × 128, or 256 × 256, and the lattice constants ax = ay = 0.1ξ. We choose the time step Δt = 0.0002 and the Josephson coupling coefficient γ = 0.3 and set the simulation temperature at T = 4.2 K. The values of the parameters coincide to the experiments in Refs. [1416]. We drop the calculation results in the early 15000 steps and start to record the effective data after that slot.

It should be known that the thermal fluctuation must be considered during simulations. The random force at each point site is independently selected from a Gaussian distribution with a zero mean. The standard deviation σd is given by[30]

where Δt is the time step and E0 is the ratio of the thermal energy to the free energy of a vortex with E0 = kBTc/ε(0)ξ(0). Here ε(0) is the free energy of a vortex per unit length at T = 0 and is defined by ε(0) = 4πξ(0)2[Hc(0)2/8π].

By using Eqs. (13)–(17) we calculate the vortex patterns at different sample sizes with different magnetic fields. Figures 2(a) and 2(b) are the vortex states we simulate at H/Hc1 = 0.4 and T = 4.2 K for the single- and two-gap superconducting systems, respectively. The Ginzburg–Landau parameter κ we select in Fig. 2(a) is κ = 3.68, which means that the system is pure type-II; while in Fig. 2(b) we choose the two parameters as κ1 = 3.68 and κ2 = 0.68 indicating that the system is type-1.5. In the following, we compare Fig. 2(a) with the vortex patterns of NbSe2 (type-I superconductor) observed by the scanning Hall probe microscopy (SHPM) at T = 4.2 K after performing an field-cooling with 1 Oe magnetic field,[16] and compare Fig. 2(b) with the vortex patterns of MgB2 (type-1.5 superconductor) observed by the scanning SQUID microscope at T = 4.2 K with 10 πT magnetic field.[15]

Fig. 2. Vortex patterns at T = 4.2K. (a) Single-gap for κ = 3.68, and (b) two-gap for κ1 = 3.68 and κ2 = 0.68 with γ = 0.3. The higher values of the Cooper pair density are shown in red.

From Figs. 2(a) we can see that the vortex lattices exhibit nearly perfectly long-range ordered, in good agreement with the experimental result of NbSe2,[16] due to the repulsive interaction between vortices. The only difference we can see is that the pattern in Fig. 2(a) is not a well-known triangular lattice as shown in the previous experiment[16] but a square one. This is because our simulation is a serious time consuming work and we cannot choose a large enough sample to overcome the confinement effect from the surface. We believe that one can obtain a triangular lattice if only he increase the sample size to a critical value. Thus, we can still speak that our simulation data are in good agreement with the experimental data in this situation.

On the other hand, for the two-gap superconducting system with κ1 = 3.68 and κ2 = 0.68, due to the fact that both the Josephson coupling condensates influence each other, vortices may attract at large separation and repel at short distance. Both the numerical [Fig. 2(b)] and the previous experimental vortex patterns[15] are different from the single-gap case. They form disordered lattices. Some local clusters and islands of magnetic flux (or vortices) and the semi-Meissner regions appear, which definitely indicate that the system is type-1.5. Our simulation also shows that if we continue increase the external magnetic field, the small vortex clusters merge with each other, to grow into bigger clusters until the whole system loses its superconductivity.

If we further increase the sample size to 25.6ξ × 25.6ξ, as it can be see in Fig. 3, besides the vortex clusters, some vortex stripes may be found in the sample. This also attributes to the competition of the short-range attraction and long-range repulsion between vortices. Here we compare the result with the vortex distribution of MgB2 (type-1.5) single crystal observed by the SHPM at T = 4.2 K after performing an field-cooling with 2 Oe magnetic field. Our numerical result and the previous experimental result[16] are in great agreement with each other. We notice that the vortex stripes in two-gap superconducting system are not straight but rather curved, i.e. they cannot be related to any crystallographic orientations of the atomic lattice. This zigzag vortex structure is quite similar to that seen in narrow superconducting ribbons with weak pinning.[38,39] Furthermore, one can also find vortices with sixfold coordination, as in an Abrikosov lattice. These features strongly suggest that the vortex clustering is not a consequence of an inhomogeneous vortex pinning. It is the nature of type-1.5 superconductors.

Fig. 3. Vortex state at H/Hc1 = 0.4 and T = 4.2 K with large sample sizes 25.6ξ × 25.6ξ for the two-gap (κ1 = 3.68 and κ2 = 0.68 with γ = 0.3) superconducting system.
4. Summary

Using the finite-element technique we theoretically model and numerically integrate the time-dependent Ginzburg–Landau equations for two-gap superconducting systems. Due to the fact that both the Josephson coupling condensates influence each other, vortices appear an interesting property that repel each other at long-range but attract at short-lange. At small magnetic field, vortex matter can be local clusters, islands, stripes, sixfold ordered lattices, or combination of them, which definitely indicate that the system is type-1.5. With further increase of magnetic field, the small local patterns merge with each other, to grow into bigger ones till the whole system loses its superconductivity. Compared to the experiments reported in Refs. [1416], to the best of our knowledge, our simulation results are the closest work so far. This provides evidence to support that the anomalous vortex distributions observed in MgB2 are in nature, which further supports the application of the type-1.5 superconductivity scenario to clean crystals of the two-gap superconductor.

Reference
[1] Tinkham M 1996 Introduction To Superconductivity 2 New York McGraw-Hill Chap. 4
[2] Huebener R P 1990 Magnetic Flux Structures of Superconductors New York Springer-Verlag
[3] Abrikosov A 1957 Sov. Phys. JETP 5 1174
[4] Ginzburg V L Landau L D 1950 Zh. Eksp. Teor. Fiz. 20 1064
[5] Kramer L 1971 Phys. Rev. 3 3821
[6] Lin S Z 2014 J. Phys.: Condens. Matter 26 493202
[7] Babaev E Carlstrøm J Silaev M Speight J M 2017 Physica 533 20
[8] Silaev M Babaev E 2011 Phys. Rev. 84 094515
[9] Diaz-Mendez R Mezzacapo F Lechner W Cinti F Babaev E Pupillo G 2017 Phys. Rev. Lett. 118 067001
[10] Babaev E 2002 Phys. Rev. Lett. 89 067001
[11] Babaev E Speight M 2005 Phys. Rev. 72 180502
[12] Carlstrom J Babaev E Speight M 2011 Phys. Rev. 83 174509
[13] Carlstrom J Garaud J Babaev E 2006 Phys. Rev. 74 134515
[14] Moshchalkov V V Menghini M Nishio T Chen Q H Silhanek A V Dao V H Chibotaru L F Zhigadlo N D Karpinski J 2009 Phys. Rev. Lett. 102 117001
[15] Nishio T Dao V H Chen Q Chibotaru L F Kadowaki K Moshchalkov V V 2010 Phys. Rev. 81 020506(R)
[16] Gutierrez J Raes B Silhanek A V Li L J Zhigadlo N D Karpinski J Tempere J Moshchalkov V V 2012 Phys. Rev. 85 094511
[17] Hicks C W Kirtley J R Lippman T M Koshnick N C Huber M E Maeno Y Yuhasz W M Maple M B Moler K A 2010 Phys. Rev. 81 214501
[18] Ray S J Gibbs A S Bending S J Curran P J Babaev E Baines C Mackenzie A P Lee S L 2014 Phys. Rev. 89 094504
[19] Kawasaki I Watanabe I Amitsuka H Kunimori K Tanida H Onuki Y 2013 J. Phys. Soc. Jpn. 82 084713
[20] Fujisawa T Yamaguchi A Motoyama G Kawakatsu D Sumiyama A Takeuchi T Settai R Anuki Y 2015 Jpn. J. Appl. Phys. 54 048001
[21] Kouchi T Yashima M Mukuda H Ishida S Eisaki H Yoshida Y Kawashima K Iyo A 2019 J. Phys. Soc. Jpn. 88 113702
[22] Huang Y Y Wang Z C Yu Y J Ni J M Li Q Cheng E J Cao G H Li S Y 2019 Phys. Rev. 99 020502(R)
[23] Mayoh D A Pearce M J Gøtze K Hillier A D Balakrishnan G Lees M R 2019 J. Phys.: Condens. Matter 31 465601
[24] Zhao Y Lian C Zeng S Dai Z Meng S Ni J 2019 Phys. Rev. 100 094516
[25] Panda K Bhattacharyya A Adroja D T Kase N Biswas P K Saha S Das T Lees M R Hillier A D 2019 Phys. Rev. 99 174513
[26] Xu C Q Li B Feng J J Jiao W H Li Y K Liu S W Zhou Y X Sankar R Zhigadlo N D Wang H B Han Z D Qian B Ye W Zhou W Shiroka T Biswas P K Xu X Shi Z X 2019 Phys. Rev. 100 134503
[27] Geurts R Milosevic M V Peeters F M 2010 Phys. Rev. 81 214514
[28] Lin S Z Hu X 2011 Phys. Rev. 84 214504
[29] Gropp W D Kaper H G Leaf G K Levine D M Palumbo M Vinokur V M 1996 J. Comput. Phys. 123 254
[30] Kato R Enomoto Y Maekawa S 1991 Phys. Rev. B 44 6916
[31] Berdiyorov G R Milosevic M V Peeters F M 2009 Phys. Rev. 80 214509
[32] Chao X H Zhu B Y Silhanek A V Moshchalkov V V 2009 Phys. Rev. 80 054506
[33] Berdiyorov G R Chao X H Peeters F M Wang H B Moshchalkov V V Zhu B Y 2012 Phys. Rev. 86 224504
[34] He A Xue C Zhou Y H 2017 Chin. Phys. 26 047403
[35] Zhitomirsky M E Dao V H 2004 Phys. Rev. 69 054508
[36] Gurevich A 2007 Physica 456 160
[37] Mazin I I Singh D J Johannes M D Du M H 2008 Phys. Rev. Lett. 101 057003
[38] Kramer R B G Ataklti G W Moshchalkov V V Silhanek A V 2010 Phys. Rev. 81 144508
[39] Karapetrov G Fedor J Iavarone M Rosenmann D Kwok W K 2005 Phys. Rev. Lett. 95 167002